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The economic and environmental benefits of laminar flow technology via reduced fuel 
burn of subsonic and supersonic aircraft cannot be realized without minimizing the 
uncertainty in drag prediction in general and transition prediction in particular. Transition 
research under NASA’s Aeronautical Sciences Project seeks to develop a validated set of 
variable fidelity prediction tools with known strengths and limitations, so as to enable 
“sufficiently” accurate transition prediction and practical transition control for future 
vehicle concepts. This paper provides a summary of selected research activities targeting the 
current gaps in high-fidelity transition prediction, specifically those related to the receptivity 
and laminar breakdown phases of crossflow induced transition in a subsonic swept-wing 
boundary layer. The results of direct numerical simulations are used to obtain an enhanced 
understanding of the laminar breakdown region as well as to validate reduced order 
prediction methods. 
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Nomenclature 

amplitude of crossflow instability mode or secondary instability mode, measured in terms of 

peak chordwise velocity perturbation 

wing chord length normal to the leading edge 

coupling coefficient for receptivity 

diameter of cylindrical roughness element 

frequency of instability oscillations (Hz) 

forcing function representing external disturbance environment 

geometry factor from receptivity theory 

gain factor associated with instability amplification 

(_ 1) 1/2 

grid index along wall-normal direction 
freestream Mach number 

N-factor (i.e., logarithmic amplification ratio) of linear crossflow instability or secondary instability 

arbitrary flow variable 

transfer function associated with receptivity 

Reynolds number based on wing chord 

time 

Cartesian velocity components aligned with (x,y,z) axes 
chordwise coordinate in direction perpendicular to leading edge 
Cartesian coordinate normal to x-y plane 

span wise coordinate, i.e., the coordinate parallel to airfoil leading edge 
chordwise wavenumber 
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spanwise wavenumber 
wing sweep angle 

spanwise wavelength of crossflow instability in millimeters in direction parallel to leading edge 

efficiency function for localized receptivity 

dummy integration variable in chordwise direction 

normalized mode shape of instability wave 

phase function 

normalized mode shape of instability wave 
radian frequency of instability wave 


Subscripts 
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init 

ins 

lb 
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crossflow direction 

initial or reference location (typically chosen to be lower branch neutral station for primary instability 

and inflow location of x/c- 0.502 for secondary instability) 

initial or reference location (used interchangeably with subscript i) 

instability wave 

lower branch neutral station for the instability mode of interest 

linear amplification phase 

nonlinear amplification phase 

roughness 

wall 

edge streamline direction 


Superscripts 
+ = wall units 


I. Introduction 

S UBSTANTIAL reductions in the fuel burn of future aircraft have been targeted both in the U.S. 1 and in Europe 2 
in order to lower the cost of air travel and alleviate the impact of aviation on the environment. A significant 
fraction of the targeted reductions in aircraft fuel burn will come from the reductions in vehicle drag. For today’s 
commercial transport aircraft, up to one half of the total drag corresponds to skin friction drag. Because laminar 
skin friction is much lower than the turbulent skin friction, flow control via delayed boundary -layer transition over 
major aerodynamic surfaces holds the potential to provide significant reductions in the overall drag. For example., 
savings in the range of 8.7-13.5 percent have been estimated to result via the implementation of laminar flow 
technology on future aircraft with a projected entry-in- service date of 2025. 3 

NASA’s Subsonic Fixed Wing (SFW) Project investigated several options to help develop robust technology 
concepts for laminar flow control and to accelerate their insertion into the operational fleet. 1 One of those concepts 
involved the delay of crossflow-induced transition over a swept aerodynamic surface by using control input in the 
form of appropriately positioned Discrete Roughness Elements (DREs) near the leading edge. 4 Most of the previous 
experimental and computational studies of the DRE concept had been carried out for low-subsonic Mach number 
configurations with modest wing-chord Reynolds numbers of up to approximately eight million, 5,6 and with pressure 
distributions that may not be optimal for wing designs for subsonic transport aircraft flying at Mach numbers 
between 0.75 to 0.90. Consequently, the studies at both Texas A & M University 7 (TAMU) and NASA 8 under the 
SFW project focused on assessing the potential capability of the DRE concept to control swept wing transition at 
transonic Mach numbers and substantially higher chord Reynolds numbers than previous applications. This work 
was continued under NASA’s Environmentally Responsible Aviation (ERA) project. 9 ' 14 Other related activities 
under the SFW project focused on in-flight characterization of crossflow receptivity to surface roughness, 15 design 
tools 16 and high Reynolds number testing for natural laminar flow (NLF) wings. 17 These activities were 
subsequently transitioned to the ERA project. The recently initiated Aeronautical Sciences (AS) project of NASA is 
aimed at developing computer-based tools and models as well as scientific knowledge that will lead to significant 
advances in our ability to understand and predict flight performance for a wide variety of air vehicles. Transition 
research under the AS project seeks to develop a validated set of variable fidelity prediction tools with known 
strengths and limitations, so as to enable sufficiently accurate transition prediction and practical transition control 
for future vehicle concepts. 18 The objective of this paper is to present the initial results obtained under this effort. 
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Figure 1. Schematic of hierarchical levels of transition prediction from a systems viewpoint. 19 The 
holistic prediction approach, based on transition as a forced response of a nonlinear system to stochastic 
forcing, appears near the top of this hierarchy. The stagewise transfer functions (i.e., gain factors) 
associated with receptivity, linear amplification, and nonlinear growth are denoted by R, G L , and G N , 
respectively. The parenthesized question mark denotes the lack of reliable methods to bridge the gap 
between physics based modeling preceding the laminar breakdown phase and some form of statistical 
turbulence modeling farther downstream. 

Given the sensitivity of swept-wing transition to wind tunnel disturbances and the cost of flight experiments, 
computational tools will need to play an important role in the assessment of DREs for subsonic transport vehicles as 
a means of risk reduction as well as to help optimize the DRE design. The DRE concept seeks to delay transition 
via controlled seeding of subdominant crossflow modes that cannot lead to transition on their own, but can keep the 
naturally dominant instability modes at bay via nonlinear modification of the basic state and/or mode competition. 
This intrinsically nonlinear control mechanism calls for a higher -fidelity prediction approach involving all of the 
relevant stages of the transition process, namely, receptivity, linear, and nonlinear growth of primary crossflow 
instabilities and secondary instability. 19,20 (See Fig. 1 for a schematic of physics based transition prediction based on 
increasing model fidelity.) The individual ingredients of this holistic approach have been available for some time 
and partial integration of those ingredients has been demonstrated in prior works. 19-26 In particular, computations of 
this type for the low Mach number configuration used during the flight experiment by TAMU at chord Reynolds 
numbers between seven and eoght million have been reported previously. 6,26,27 Results based on a combination of 
nonlinear parabolized stability equations (NPSE) and linear secondary instability theory (LSIT) showed that the 
most unstable stationary crossflow mode is weakened in the presence of artificially introduced roughness elements, 
with an accompanying reduction in the growth potential of the high-frequency secondary instability modes that are 
likely to cause the onset of laminar-turbulent transition. 26 However, the same results also indicated that the control 
mode amplitudes required to stabilize the dominant modes can become excessively large, especially when the 
maximum N-factors of the primary crossflow instability are rather large, so that the possibility of premature 
transition due to overcontrol cannot be ruled out. Similar findings had previously been obtained for a supersonic 
swept wing configuration at Mach 2.4 and chord Reynolds number of 16 million. 23,24 

The aforementioned SFW and ERA efforts targeting passive laminar flow over a subsonic transport involved a 
series of conceptual design studies of a laminar flow control (LFC) flight test. The initial studies under SFW 7,8 
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focused on conditions relevant to a subsonic transport, namely, a flight Mach number of 0.75, a chord Reynolds 
number, Re c , between 15 and 20 million, and a leading-edge sweep angle, 77 of 30°. The resulting wing design was 
predicted to have a significant potential for natural laminar flow over both suction and pressure surfaces at the 
design angle of attack of zero degrees. At off-design conditions, however, earlier transition was deemed likely as a 
result of stationary crossflow modes. This finding opened up the possibility of DRE-based transition control on the 
same airfoil, either to increase the lift coefficient or to extend the angle -of-attack range for laminar flow. 7 The 
companion study by Li et al. 8 sought to assess the feasibility of achieving DRE-based transition control by using 
NASA’s prediction tools for nonlinear development and secondary instability stages of the transition process and to 
help mature those prediction tools and the methodology to apply them in a practical context. The latter results are 
extended in this paper (as well as in a companion paper 27 ) to address the receptivity phase (i.e., roughness -induced 
excitation of stationary crossflow modes), non-parallel and nonlinear effects on the growth of high-frequency 
secondary instabilities of stationary crossflow modes, and the subsequent breakdown process. 

In the interest of making this paper self-contained, a brief description of the flow configuration as well as the 
computational methodology is given in Section II below. Results pertaining to the excitation of stationary crossflow 
vortices due to deterministic but spatially distributed surface roughness (waviness) are outlined in Section III. 
Section IV presents the findings from direct numerical simulations aimed at validating the previous predictions 
based on both LSIT and parabolized stability equations (PSE) which pertain to the amplification of secondary 
instability of stationary crossflow vortices. Section V outlines a specific scenario for the nonlinear breakdown of 
secondary instability modes resulting in a turbulent, swept- wing boundary layer. To gain a further understanding of 
the role of traveling crossflow modes in the transition of a swept wing boundary layer, the secondary instability 
characteristics of finite amplitude traveling crossflow modes are described in Section VI. Concluding remarks are 
presented in Section VII. 


II. Flow Configuration and Analysis Codes 

The design process for the laminar flow airfoil TAMU-003T-75(v.90) was described by Belisle et al. 7 The 9.3 
percent thick, 30-degree swept airfoil is designed to achieve natural laminar flow over 60 and 50 percent of the 
suction and pressure surfaces, respectively, at the design condition of zero degree angle of attack, Mach number of 
M = 0.75, and chord Reynolds number of approximately 17 million based on freestream speed, a streamwise chord 
length of 3.658 m (12 ft), freestream temperature of 216.65 K and density of 0.302 kg/m 3 . The latter conditions 
correspond to a flight altitude of 12,192 m (40,000 ft). Design constraints for the wing included a lift coefficient 
that is typical of subsonic transports and a wing thickness distribution that is suitable for a mid- size business jet. At 
the design angle of attack, both Tollmien-Schlichting and crossflow instabilities are sufficiently weak that natural 
laminar flow should be achievable over a significant portion of the airfoil surface (x/c>0.6 on the suction surface and 
x/c>0.5 on the pressure side) without any external means of boundary layer control. In this paper, we focus on the 
off design condition corresponding to an angle of incidence equal to -1 degree, i.e., the same condition as that used 
by Li et al. 8 in their PSE and LSIT analysis. At this off-design condition, a stronger crossflow instability exists on 
the suction side, 7 ’ 8 which provides the basis for potentially implementing DRE-based laminar flow technology at the 
higher Reynolds numbers of interest. 

For the work described in this paper, the mean boundary -layer flow over the suction surface of the airfoil was 
obtained from the boundary layer solution obtained by Li et al. 8 It was computed using a boundary layer solver 
BLSTA, 28 designed specifically for generating accurate boundary layer solutions that are required for boundary layer 
stability analysis of swept and tapered wings, in conjunction with the inviscid surface-pressure distribution derived 
from an Euler solution under free flight conditions. 7 The number of wall-normal points used for boundary layer 
computations was well in excess of that used by Wie 28 for validating the BLSTA code. Sensitivity to streamwise 
resolution was also assessed by Li et al. 8 by using 115, 229, 457 and 913 points over the chord length of interest 
(with proportionally greater number of points near the leading edge), and the solution was converged with 115 
streamwise points. Similar to the preceding work, 8 the computations presented in this paper are based on 229 points 
over the chord length. The boundary layer solution computed in this manner is suitably interpolated using cubic 
splines when necessary. 

The growth of high-frequency secondary instability modes supported by the finite amplitude stationary crossflow 
vortex was analyzed by Li et al. 8 with a fully spatial methodology described in earlier papers 23,24 and their 
predictions for the case of Ainit = lxlO' 5 are compared with the DNS results in Section IV. The main difference 
between the primary and secondary stability analyses is that the basic state for the secondary modes (i.e., the mean 
boundary layer flow modified by the primary crossflow mode) varies in both surface normal and spanwise 
directions; and hence, the instability characteristics of the secondary modes must be analyzed using a planar, partial 
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differential equation based eigenvalue problem, rather than an ordinary differential equation based eigenvalue 
problem of the classical analysis. The selection of grid and other aspects of the numerical solution was based on 
extensive experience with a similar class of flows and spot checks were made to ensure that the impact of reasonable 
variations with respect to those choices was negligible. Typically, 121 points in the wall-normal direction and 33 
points in the span wise direction were used for the 2D eigenvalue analysis. A similar methodology and spatial 
resolution are used in Section VI to analyze the secondary instability of traveling crossflow modes. 

Linear and nonlinear development of the primary instability (i.e., either stationary or traveling crossflow modes) 
was computed using linear/nonlinear parabolized stability equations (LPSE and NPSE, respectively) as implemented 
in the Langley Stability and Transition Analysis Codes (LASTRAC). 29 For nonlinear PSE marching, 281 points are 
used in the wall-normal direction (which is more than sufficient), and the stream wise marching step-size was 
appropriately adjusted for each computation, if necessary, to capture the corresponding modal evolution. We denote 
the ranges of Fourier modes in the span wise direction and in time by (-m, m) and (-n, n ), respectively. A total of 48 
span wise Fourier harmonics (m = 48) is used to compute the nonlinear development of a primary crossflow mode of 
specified spanwise wavelength, A, and the truncation error is estimated a posteriori to be at least 9 orders of 
magnitude smaller than the most energetic fundamental mode. 

Receptivity calculations are performed using the adjoint PSE module implemented within the LASTRAC code. 30 
The spatial resolution used for the adjoint PSE computations is similar to that used for the regular LPSE 
computations. 

To study the breakdown of finite-amplitude crossflow vortices due to high-frequency secondary instability, a 
direct numerical simulation (DNS) is performed by solving the Navier-Stokes equations in generalized curvilinear 
coordinates. The working fluid is assumed to be an ideal gas with linear Newtonian strain-stress relation. The 
Fourier law is used for the heat flux. The computational domain is approximately aligned with the axis of the 
primary crossflow mode and extends over more than 100 wavelengths of the secondary instability wave at a 
frequency of 42.5 kHz. The size of the computational grid corresponds to 5,000x242x307 points in the 
chordwise, spanwise, and wall-normal directions, respectively. In terms of wall units of the fully turbulent flow near 
the outflow of the domain, the grid resolution in the streamwise (v) and spanwise directions (y) corresponds to Av + « 
10 and Ay + « 4.5, respectively. The first grid point in the wall normal direction is located at A z + ~ 0.6. 

A 7 th -order weighted essentially non-oscillatory (WENO) scheme 32,33 is used to compute the convective flux 
terms. This particular WENO scheme combines a high order of accuracy with relatively low dissipation, making it 
suitable for simulations of compressible transitional flows. The resolution properties of this scheme are documented 
in many references, e.g. Martin et al. 34 For the viscous flux terms, a 4 th -order central difference scheme is used and 
the 3 rd -order low storage Runge-Kutta scheme by Williamson 35 is employed for time integration. Unsteady non- 
reflecting boundary conditions based on Thompson 36 are imposed at the freestream boundary. At the outflow 
boundary, a buffer domain technique is used before the primitive variables are extrapolated at the outflow plane. The 
DNS code has been previously shown to be suitable for computing fully turbulent flows. 34,37 In additional 
unpublished work, further validation of this code was performed in the context of the propagation of linear 
instability waves in a 2D, high-speed boundary layer. Convergence studies showed that the 7 th -order WENO 
performs well at a resolution of 10 points per wavelength or higher, depending on the streamwise extent of the 
computational domain. 


III. Receptivity Characteristics 

Even though the N-factors for stationary crossflow vortices are lower than those of traveling crossflow modes, the 
inevitable imperfections in surface finish lead to preferential excitation of stationary modes. 38,39 As a result of their 
higher initial amplitudes, these stationary modes tend to dominate the crossflow transition process in a benign 
disturbance environment analogous to that in flight. For cases where the roughness amplitude is sufficiently small 
and the wavenumber spectrum of the surface height perturbations overlaps with the length scales of the crossflow 
vortex instabilities, the resulting receptivity can easily be predicted by extending the Goldstein-Ruban theory 40,41 for 
localized receptivity. Thus, the flow perturbations associated with an arbitrary crossflow mode of spanwise 
wavenumber, /?, can be expressed in the form 42 : 

q(x, y, z\ o), p) = C(x,a>,fi)(/)(x,y) exp(/(^ ad% + pz-cot) (1) 

where the cumulative coupling coefficient for distributed receptivity, C, corresponds to the effective instability 
amplitude at a selected reference location v = xq. The quantity $x,y) denotes a normalized mode shape associated 
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with the perturbations produced by the instability mode. As a result of the assumed linearity of the receptivity 
process, the coupling coefficient can be expressed as: 

C= VF *o ( ’ ' )exp( iins{ ' ))d (2) 

ins (x ' > = x 0 ins < ’ > d W 

where A(£, P) characterizes the efficiency of the local receptivity to surface roughness, i.e., the local Green’s 
function with the response measured at the forcing location, The geometry factor F represents the Fourier 
coefficient of the surface roughness distribution corresponding to the span wise wavenumber /?. The decoupling 
between geometry and efficiency functions was first pointed out by Goldstein in his pioneering work for localized 
excitation of Tollmien-Schlichting waves. 40 Computing the efficiency function for a given mean flow and wave 
number thus provides a basis for further predicting the disturbance amplitude associated with a surface non- 
uniformity of finite spatial extent. 



Fig. 2. Localized receptivity due to a single array of roughness dots mounted near the leading edge of the 
TAMU-003T-75(v.90) airfoil: Effect of array location, x n on the amplitude of stationary crossflow vortex. 

(A = 2.67 mm; D r = 1.33 mm; and the height of roughness dots is assumed to be 1.0 pm.) 

For the TAMU-003T-75(v.90) airfoil configuration, 7 the stationary crossflow vortices that achieve the highest 
linear amplification at the off design condition of interest correspond to a spanwise wavelength of approximately 
eight mm. 7,8 The A = 2.67 mm mode has a linear N factor of approximately six and was shown to be suitable as a 
control input based on DREs. 8 The excitation of the control mode using a spanwise periodic array of cylindrical 
roughness elements (i.e., DREs) of diameter D r = 1.33 mm is considered first. Fig. 2 shows the magnitude of the 
coupling coefficient plotted as a function of the chord wise location v r of the array of roughness elements. The peak 
effective receptivity occurs when the DREs are placed near x/c = 0.0033, i.e., near the lower branch location for A = 
2.67 mm. 

The stationary crossflow mode with A = 8 mm is excited as part of the flow response to naturally occurring 
surface roughness, which in principle may include both deterministic and random components. As described by 
Choudhari, 43 the receptivity to random distributed roughness can be evaluated using the solution for receptivity to a 
wavy wall as a building block component. Computational predictions for receptivity to naturally occurring, 
stochastic roughness distribution over a flight article were presented by Carpenter et al. 44 Receptivity to wavy walls 
was also considered in the latter paper, but those results were restricted to the canonical problem of a swept 
Hiemenz boundary layer. As previously shown by Collis and Lele, 45 the receptivity to localized roughness elements 
is strongly influenced by the curvature of the airfoil surface. Thus, it is useful to examine the distributed receptivity 
due to wavy undulations of small amplitude in the airfoil surface. Figure 3 shows the coupling coefficient for a 
stationary crossflow vortex with A = 8 mm, plotted as a function of the chordwise wavenumber, a w , of the wavy 
surface height distribution. Again, similar to the canonical case of a swept Hiemenz boundary layer, the peak 
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receptivity levels are confined to a narrow band of surface wavenumbers. The wavenumber corresponding to the 
maximum coupling coefficient is close to the neutral wavenumber of the stationary crossflow mode with 2 = 8 mm, 
but is slightly detuned from it. Similar to the swept Hiemenz flow, the wavenumber bandwidth of the flow response 
is substantially larger than that of acoustic excitation of T-S waves through scattering by the surface waviness. 42 



Fig. 3. Receptivity of boundary layer flow over the 
TAMU-003T-75(v.90) swept airfoil configuration to 
wavy wall roughness height distribution with peak 
height of 1.0 pm: Effect of surface wavelength on the 
amplitude of stationary crossflow vortex with 2 = 8.0 
mm. 



Fig. 4. Distributed receptivity due to wavy wall 
compared with receptivity due to a single array of 
roughness dots at the lower branch neutral station 
(D r = 1.33 mm and the height of the dots is assumed 
to be the same as the amplitude of the wavy wall). 


To ascertain the relative receptivity levels for surface waviness relative to those for a single isolated region of 
surface imperfection, the cumulative coupling coefficient for a w = a w> max = 363.2 m" 1 (which corresponds to the 
peak of the response curve in Fig. 3) is plotted as a function of chordwise location in Fig. 4. For reference, the 
coupling coefficient for a local imperfection corresponding to a single array of cylindrical roughness elements 
centered on the lower branch station is also shown in this figure. The diameter of the dots is assumed to equal one 
half of the spanwise wavelength. For a w = a w> max , the stationary crossflow perturbations due to portions of the 
wavy surface in the vicinity of the lower branch tend to reinforce each other with a nearly zero cancellation. Thus, 
the cumulative response is approximately one order of magnitude larger than that of a single array of roughness 
elements centered at the neutral station. The resulting values of coupling coefficient indicated that the initial 
crossflow amplitudes considered in the NPSE calculations of Li et al. can easily be excited by a wavy wall height 
distribution at the peak wavenumber even if the root-mean-square height of this distribution is just a few microns. 

IV. Non-Parallel and Nonlinear Effects on Secondary Instability Amplification 

Li et al. 8 performed NPSE computations to predict both isolated and combined nonlinear evolution of the target 
(2 = 8 mm) and control (2 = 2.667 mm) modes of the primary, i.e., stationary crossflow instability for varying initial 
amplitudes. To obtain a greater insight into the expected transition locations, secondary instability analysis was 
performed using the NPSE results as the modified basic state. Both LSIT and PSE analyses were used to 
characterize the N-factor evolution of secondary instability modes over a broad range of frequencies. In most cases, 
the modified basic state supported multiple families of secondary instabilities that derived their amplification from 
turbulence production mechanisms associated with wall-normal and/or spanwise shear of the basic state. In the 
notation of Ref. 8, these modes are referred to as the Y- and Z-modes, respectively; and we will continue to employ 
that nomenclature in this paper. For 2=8 mm, the secondary instability modes associated with wall-normal shear 
(i.e., Y-mode in the notation of Ref. 8) had the higher amplification rates and, therefore, might be expected to 
dominate the laminar breakdown process over the swept airfoil configuration examined in this paper. This 
breakdown process is examined in a companion paper by Duan et al. 27 

The receptivity mechanisms associated with the secondary instability modes are not well understood at this 
point; and therefore, it is conceivable that, under certain circumstances, the secondary instability modes with lower 
amplification factors might get excited more strongly than the modes with dominant linear growth and hence cause 
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the onset of transition. Indeed, there is some evidence based on prior experiments 46 and computations 47 that the 
secondary instability modes associated with span wise shear (i.e., Z-mode) may become relevant to the transition 
process even when their amplification factors are lower than the wall -normal- shear modes. Furthermore, for the 
airfoil configuration and disturbance parameters examined by Malik et al., 14 the Z-mode of secondary instability had 
higher peak amplification factors than the Y -mode. Thus, it is also important to consider the linear and nonlinear 
evolution of the Z-modes. 

According to the LSIT analysis of Li et al., 8 a span wise shear mode of frequency 42.5 kHz achieved peak linear 
N-fators approaching seven and, therefore, could cause transition if it had a sufficiently large initial amplitude. To 
validate the LSIT and PSE predictions, DNS computations were performed for two different initial amplitudes that 
differed by a factor of more than 150. The smaller initial amplitude (A init = 3.4xl0‘ 5 ) was chosen to mimic the linear 
growth of the secondary instability mode, whereas the higher initial amplitude (A init = 0.0055) was chosen to capture 
nonlinear modification of the initial growth, before a breakdown to turbulence. Fig. 5 displays the amplitude 
evolution (based on the peak chordwise velocity fluctuation over the cross-section of the primary crossflow vortex) 
from the two DNS computations and their comparison with the LSIT and PSE predictions of Li et al. 7 The DNS 
result for the smaller initial amplitude agrees very well with the PSE prediction, confirming that the differences 
between the LSIT and PSE results are caused by neglecting the weakly nonparallel effects in the LSIT predictions. 
Furthermore, the DNS results for the higher initial amplitude confirm the initially stabilizing effect of nonlinearity 
on the fundamental mode amplitude, which is analogous to the finding of Li et al. 8 The approach to quasi-saturated 
behavior in the fundamental mode evolution is followed by a rapid increase in the modal amplitude and then a 
stochastic behavior that is suggestive of the onset of transition. The details of this transition zone are further 
explored in Section V below. 



X/C 

Fig. 5. Amplification curve for Z-mode of secondary instability: Comparison of DNS solutions with LSIT 
and PSE predictions 

V. Laminar Breakdown due to Secondary Instability of Stationary Crossflow Vortices 

Li et al. 8 performed NPSE computations to track the joint nonlinear evolution of a stationary crossflow mode 
with /l = 8 mm and A init (xi b ) = lxl0‘ 5 and its secondary instability modes of varying initial amplitudes upto the 
location where the mean skin friction began to rise rapidly. However, the details of the transition zone could not be 
inferred due to the limitations of the PSE approach. In this section, the results of the DNS computation involving 
the same primary mode and a span wise-shear type (i.e., Z-mode) secondary instability with /= 42.5 kHz and A init 
( x/c = 0.502) = 0.0055 are used to explore the boundary layer evolution through the laminar breakdown region. Due 
to the large computational grid and correspondingly long time interval required to average over a large number of 
eddies, fully converged turbulence statistics have not been acquired thus far. However, a comparison of 
instantaneous, spanwise-averaged statistics at different instants of time indicates that the primary conclusions 
presented herein may be expected to remain unchanged with further accumulation of time history data. Unless 
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otherwise specified, all of the flow statistics presented herein have been averaged over 20 periods of the 
fundamental secondary mode. 



Fig. 6. Dimensional mean wall velocity at first point off the wall (which is located at a fixed wall-normal 
distance over the entire length of the computational domain). Due to averaging of flow statistics over a 
limited time interval, additional moving averaging in x has been performed to smooth short scale 
oscillations in the evolution of flow statistics. Subscripts rand c indicate velocity components aligned with 
the edge streamline and crossflow direction, respectively. 



Fig. 7. Instantaneous chordwise velocity contours 
over both the horizontal coordinate plane k = 3 
(which lies very close to the wall) and the vertical 
plane x/c = 0.537 (which lies just upstream of the 
onset of transition). 



Fig. 8. RMS fluctuations in chordwise velocity 
component across the cross-section of the 
stationary vortex at selected stations (x/c <0.544). 
The RMS is taken over 20 periods of the 
secondary instability wave. 


The chordwise evolution of near-wall mean velocities at a fixed distance above the surface (Fig. 6) indicates a 
rapid rise in the near-wall mean velocity close to x/c = 0.54, indicating the onset of transition location. The skin 
friction reaches a peak near x/c = 0.55 and then decreases weakly over the remaining length of the computational 
domain. No prominent overshoot in skin friction is observed in this particular case. Whether the absence of an 
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overshoot in the skin friction is caused by the spanwise averaging to obtain the result in Fig. 6 or whether it is a 
genuine feature of the present scenario for stationary crossflow breakdown remains to be seen. Flow statistics over a 
longer time interval will be used in the future to clarify whether the overshoot is also absent from the chordwise 
evolution of mean wall shear within individual spanwise planes. Fig. 6 indicates that the wall shear associated with 
the crossflow profile also increases after the onset of transition. However, the magnitude of this increase is weaker 
in comparison with the wall shear along the direction of the inviscid streamline. Moreover, the eventual decrease in 
the wall shear associated with the crossflow profile implies that the wall shear direction approaches the direction of 
the streamline at the edge of the boundary layer as the flow goes through transition. The latter finding, originally 
noted by Huai et al. 48 appears to be associated with the decrease in both crossflow magnitude and the crossflow 
shear, which corresponds to a shift towards increased two -dimensionality in the mean turbulent flow. 

Fig. 7 displays the spanwise variation in instantaneous contours of chordwise velocity distribution over the k=3 
surface, which is located at a fixed distance from the airfoil surface. The steep gradients in chordwise velocity 
contours resemble the sawtooth transition front that is commonly seen in naphthalene -based surface flow 
visualizations of crossflow transition dominated by stationary vortices. The origin of the turbulent wedge lies within 
a spanwise region that corresponds to relatively smaller boundary layer thickness, and the completion of the 
transition process (i.e., the location of peak mean wall shear in Fig. 6) approximately coincides with the merging of 
adjacent pairs of turbulent wedges. The half-angle of the turbulent wedge inferred from the spanwise wavelength of 
the stationary crossflow vortex and the cord wise extent of breakdown region corresponds to approximately six to 
seven degrees. A similar sawtooth pattern has also been observed by Duan et al. 27 in their DNS of transition caused 
by the Y-mode of secondary instability (namely, the mode that is primarily associated with the wall -normal shear of 
the basic state). The similarity between the sawtooth patterns observed in both simulations indicates that the 
mechanism behind the breakdown of stationary crossflow vortices cannot be inferred from the shapes of the 
respective transition fronts. 

The evolution of (short duration) root mean square (RMS) chordwise velocity contours across the cross-section 
of the stationary crossflow vortex is shown in Figs. 8 through 10. The initial evolution in Fig. 8 is qualitatively 
analogous to that presented by Li et al. 8 on the basis of their NPSE computations. When transition is completed (see 
image 4 in Fig. 9), the peak fluctuations have shifted to the nearwall region. However, a visible signature of the 
stationary vortex is still observed within the outer part of the boundary layer. The spanwise variations appear to 
become more diffuse as the turbulent flow evolves further downstream; yet they can still be observed even at x/c = 
0.578 (Fig. 10), i.e., more than 10 boundary -layer thickneses downstream of the peak mean wall shear. 



0.556 


Fig. 9. RMS fluctuations in chordwise velocity 
component across the cross-section of the 
stationary vortex at selected stations (0.544 < x/c < 
0.556). 



Fig. 10. RMS fluctuations in chordwise velocity 
component across the cross-section of the stationary 
vortex at selected stations (0.556 < x/c < 0.578). 


Figure 1 1 shows an instantaneous snapshot of the flowfield at a selected instant of time. An isosurface of u = 
175 m/s (i.e., approximately 71 percent of the edge velocity) is shown along with crossplane visualizations of 
temperature contours at three selected chordwise stations. Upstream of the transition onset, the isosurface indicates 
a series of ridges that are inclined with respect to the vortex axis and are centered along a spanwise location that 
approximately corresponds to the peak location of the secondary instability mode shape. The latter location lies 
approximately in-between the adjacent pair of locations where the boundary layer thickness has increased 
significantly as a result of the upwelling associated with a finite -amplitude crossflow vortex. As nonlinear effects 
come into play, the ridges develop a localised spanwise oscillation near the peak of the secondary instability mode 
shape (see the middle cut at x/c = 0.541 in Fig. 11). Further downstream, this local deformation begins to spread 
out, particularly towards the wall and at approximately the same spanwise location relative to the axis of the 
stationary crossflow vortex. Since the peak of the Z-mode of secondary instability was below the crests of the 
overturning contours, the tertiary fluctuations can quickly begin to influence the nearwall dynamics, which then 
leads to the onset of self-sustaining turbulence. 
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Fig. 11. Instantaneous isosurface of u - 175 m/s along with temperature contours at selected x 
locations (. x/c = 0.527, 0.541, and 0.549, respectively). The chordwise velocity for the isosurface is 
selected to coincide with the approximate location of peak chordwise velocity fluctuations associated 
with the secondary instability mode of interest. 


VI. Secondary Instability of Traveling Crossflow Vortices 

In general, a swept-wing boundary layer can be susceptible to various types of primary instabilities such as 
attachment-line instability, 49 stationary and traveling crossflow modes, and Tollmien-Schlichting (TS) waves. 
Regardless of which mode dominates the primary amplification stage, the other modes could play a role during the 
nonlinear stage and, hence, influence the onset of transition. Although it appears that the more likely scenario in 
cases of practical interest would involve the breakdown of stationary crossflow vortices, either via some type of high 
frequency secondary instability 20,50 ’ 53 or through interactions involving TS waves or the traveling mode of crossflow 
instability, other scenarios involving a more significant role for the traveling crossflow instability 54,55 may become 
relevant under less common circumstances. These may happen when either the initial amplitudes of traveling modes 
become sufficiently large as a result of higher amplitude freestream disturbances and/or when the surface has been 
carefully polished to minimize the initial amplitudes of stationary crossflow modes. 

In a practical setting, it is perhaps more important to investigate the nonlinear breakdown scenarios associated 
with a mixed evolution of stationary and traveling crossflow modes. However, as a prelude to the more complex 
scenario, the secondary instability of traveling crossflow modes alone is studied herein to obtain some useful 
insights regarding their breakdown. 

A. Nonlinear Evolution of a Traveling Crossflow Mode 

Figure 12 shows N-factor curves of crossflow instability of 8 mm span wise wavelength for selected frequencies 
ranging from 0 to 2,250 Hz. As expected, the most amplified traveling mode achieves higher N-factors than the 
stationary mode; and, more importantly, the traveling mode with /= 1,500 Hz reaches N = 10 within nearly half the 
distance it takes for the stationary crossflow modes to reach the same N-factor value. This indicates the utility of 
studying the breakdown of traveling waves to obtain a more complete picture of crossflow-dominated transition, 
especially in noisy environments. Since an N-factor value of 10 is first reached by a traveling crossflow mode of f- 
1,500 Hz, the latter mode is selected for the purpose of secondary instability analysis. 
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Figure 12. Linear N-factors for crossflow instability modes with A, = 8 mm and frequency ranging from 0 to 
2,250 Hz. The N-factor first reaches a value of 10 for the 1,500 Hz mode indicated by green curve. N -factor 
curve for the stationary crossflow mode is shown in red. 

Figure 13(a) compares the nonlinear development of stationary and traveling crossflow vortices of the same 
span wise wavelength (k = 8 mm) and the same initial amplitude (A init = lxlO' 5 ). In the initial region of linear 
growth, the amplitudes of modes associated with the traveling crossflow vortex (green lines in Figure 13(a)), rise 
much faster than the stationary mode amplitudes (red lines in Figure 13(a)). As the amplitudes become larger, the 
traveling mode begins to saturate much earlier than the stationary mode, possibly due to the rapid rise in its mean 
flow correction, and ceases its growth beyond approximately x/c = 0.2. On the other hand, the stationary mode does 
not saturate until approximately x/c = 0.5. The traveling mode reaches the peak amplitude of approximately 0.16 
before beginning to decay in an oscillatory fashion. On the other hand, the stationary mode grows monotonically 
over most of the stream wise domain and reaches a maximum amplitude of 0.32, i.e., twice the peak amplitude 
achieved by the traveling mode. The above comparison between nonlinear evolution of isolated traveling and 
stationary crossflow disturbances is very similar to that reported by Malik et al. 14 for a different airfoil configuration. 
The maximum mean flow correction amplitudes for both stationary and traveling mode are relatively similar at 0.08 
and 0.07, respectively. 

Modal amplitudes of traveling crossflow vortices of 1500 Hz are shown in Figure 13(b) for initial amplitudes of 
2.5x1 O' 6 , 5. 0x1 0' 6 , 7.5x1 O' 6 and lxlO' 5 . Amplitudes of the mean flow correction and the first harmonic are also 
shown for the case of the smallest initial amplitude. The general features of oscillatory behavior and decay after 
saturation as seen earlier in Figure 13(a) for A init = lxlO 5 are present for all initial amplitudes, with waves 
corresponding to higher initial amplitudes exhibiting an earlier quasi-saturation. Figure 13(b) also indicates that 
traveling crossflow modes with higher initial amplitudes saturate at progressively smaller peak amplitudes. 

When the initial amplitudes of the traveling and stationary crossflow amplitudes are equal, the stronger growth 
of the traveling mode leads to an earlier onset of secondary instability. Figure 14 presents a comparison of the 
secondary N-factors for the stationary and traveling primary modes considered in Figure 13(a). Specifically, the N- 
factor curves of Y- and Z-modes of secondary instability that first reach an N-factor of 10 are plotted for both types 
of primary waves. The peak N-factor for the Z-mode instabilities of the stationary crossflow vortex is less than 10; 
and, hence, the secondary Z-mode disturbance with f- 42.5 kHz corresponding to the peak of the N-factor envelope 
is shown in that case. In all cases the neutral points of the secondary instability appear to be in the neighborhood of 
the locations where the primary crossflow vortices begin to saturate (see Figure 13). The most dangerous secondary 
instability mode of the traveling crossflow vortex (i.e., a Y-mode) becomes unstable at a primary vortex amplitude 
of approximately 0.16. Similarly, the corresponding mode in the stationary case (also a Y-mode) becomes unstable 
at nearly the same primary amplitude of approximately 0.17. Over the same distance of approximately 20% chord 
length, the most unstable Y-modes reach an N-factor of over 30 in both the stationary and the traveling cases. 
However, the initial slope of the Y-mode N-factor curve for the traveling vortex is larger than that in the stationary 
case, therefore, the secondary instability reaches an N-factor of 10 faster than than in the stationary case. 

In the future, it will be useful to extend the secondary instability computations of traveling crossflow modes to 
even lower initial amplitudes, so as to provide insights regarding the threshold initial amplitude of traveling 
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crossflow modes that could dominate the transition process. 



Fig. 13(a). Nonlinear evolution of stationary and 
traveling crossflow (f = 1500 Hz) modes with X = 8 
mm and the same initial amplitude (A init = lxlO' 5 ): 
Comparison of modal amplitudes for the 
fundamental and mean flow correction modes. 

B. Secondary Instability of traveling crossflow vortex 



Figure 14. Comparison of N-factor curves for most 
dangerous secondary instability modes of 
stationary and traveling crossflow vortices. The 
metric for growth potential corresponds to earliest 
attainment of N = 10 or the highest N-factor when 
N = 10 is never achieved. Both primary crossflow 
vortices have the same spanwise wavelength of 8 
mm, and the traveling crossflow vortex has a 
frequency of 1500 Hz. 



Fig. 13(b). Evolution of fundamental amplitude 
for traveling crossflow vortices with different initial 
amplitudes. The amplitudes of mean flow correction 
and the first harmonic are also shown for the case 
with the smallest initial amplitude (A init = 2.5xl0' 6 ). 



Fig. 15. Peak N-factors of secondary instability of 
traveling crossflow vortices of different initial 
amplitudes. 


Figure 15 shows the N-factor curves for secondary instability modes that first achieve an N-factor of 10 for each 
initial amplitude of the traveling crossflow vortex, namely A init = 2.5xl0' 6 , 5.0xl0' 6 , 7.5xl0' 6 , and lxlO' 5 , 
respectively. Regardless of the initial primary amplitude in this case, the dominant secondary instability mode is a 
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Z (mm) 

Fig. 16(a). Chordwise velocity eigenfunction for a Y-mode secondary instability of a primary traveling 
crossflow wave of initial amplitude lxlO' 5 . Secondary instability frequency is 67.5 kHz. 



Z (mm) 

Fig. 16(b). Chordwise velocity eigenfunction for a Y/Z mode secondary instability of a primary traveling 
crossflow wave of initial amplitude lxlO' 5 . The frequency of the secondary instability mode is 20 kHz. 



(c) A init = 7.5xl0' 6 (d) A inU = l.OxlO' 5 


Fig. 17. Chordwise velocity contours associated with a finite amplitude traveling crossflow mode. 
The coordinate system is assumed to translate with the traveling mode and the chordwise location 
of x/c = 0.20 is chosen to be in the vicinity of the onset of secondary instabilities. 

Y-Mode except for the case of lowest initial primary amplitude for which no Y-mode was found by the eigenvalue 
solver. The subdominant modes are of the mixed type, i.e., Y/Z modes that receive comparable contributions from 
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both wall-normal and spanwise shear of the basic state. For the cases with larger values of initial primary amplitude, 
the peak N-factors values for secondary instability are correspondingly larger as expected. While the Y -modes decay 
rapidly after reaching their peak N-factors, the Z-modes continue to grow almost monotonically. The frequency of 
Y-mode disturbances achiving the highest N-factor values ranges from 58-67.5 kHz, whereas the frequency range 
for the most amplified Y/Z modes corresponds to 15-19 kHz. 

Figs. 16(a) and 16(b) illustrate the chordwise velocity eigenfunctions of the Y mode and Y/Z mode, respectively, 
at x/c = 0.25. In a coordinate frame that travels with the primary instability (i.e, the traveling crossflow mode), the 
basic state becomes steady in time. Figure 17(a) through 17(d) show the mean chordwise velocity contours in the 
traveling frame at x/c = 0.2 for initial primary amplitudes of 2.5xl0' 6 , 5.0xl0' 6, 7.5xl0' 6 and 10' 5 , respectively. For 
the highest initial amplitude, the base velocity contours develop an overturning behavior similar to that seen in 
previous studies of finite amplitude stationary crossflow vortices. Given that the onset of secondary instability even 
for the lowest initial amplitude (2.5x1 O' 6 ) occurs near x/c = 0.2, it may be inferred that even a relatively modest 
deformation of the underlying basic state can lead to the amplification of secondary modes. 

VII. Concluding Remarks 

In this paper, previously reported holistic predictions of crossflow induced transition over a realistic swept airfoil 
configuration were extended via high fidelity computations involving the receptivity, secondary instability, and 
nonlinear breakdown stages of transition. The results presented herein demonstrate that the adjoint PSE capability 
of the LASTRAC code can be used to predict the excitation of stationary crossflow modes due to small amplitude 
surface roughness for a realistic airfoil configuration of relevance to future subsonic transports. Furthermore, the 
compressed length of the stable (i.e., subcritical) region upstream of the lower branch neutral location near the 
airfoil nose is shown to reduce the peak receptivity associated with nonlocalized, wavy-wall distributions of surface 
height variation. In future, the capability for receptivity prediction will be further extended to allow more rational 
predictions of the overall transition process for a given specification of the environmental disturbances. 

Despite the substantial differences in mode shapes of the secondary instability modes associated with spanwise 
and wall-normal shear, the overall features of the transition pattern obtained herein for the spanwise mode are very 
similar to those obtained in a separate work focused on transition due to secondary instability modes driven by the 
wall-normal shear. Specifically, the transition front visualized via time averaged wall shear indicates a spanwise 
varying pattern resembling the sawtooth pattern visualized in previous wind tunnel experiments using naphthalene 
sublimation. Furthermore, in both cases, the memory of the initial stationary crossflow vortex persists even after 
transition. In the present computational domain, distinct spanwise inhomogeneity is observed even at distances of 
more than 20 boundary layer thicknesses past the location of peak skin friction. In other words, the ensemble- 
averaged turbulent flow varies along the spanwise direction throughout the computational domain. Since a drag 
prediction method based on the assumption of spanwise homogeneity of the first and second moment statistics will 
not capture these 3D effects, skin friction predictions based on Reynolds Averaged Navier-Stokes (RANS) equations 
may entail a nonzero bias in its predictions. Quantification of these effects associated with the relaxation of the 
spanwise perturbed turbulence will be a topic for future research. 

The traveling crossflow vortices support secondary instabilities in much the same way as the stationary 
crossflow modes. Even though the traveling crossflow vortices saturate at much lower amplitudes than their 
stationary counterparts, the growth of secondary instabilities of the traveling crossflow vortices is at least as strong 
as the secondary modes of the stationary crossflow modes. The important issue of threshold initial amplitudes of the 
traveling vortices required to sustain a significant amplification of secondary instabilies will also be examined in 
future work. 

The multiplicity of secondary instability modes with a significant growth potential implies that transition 
prediction based on secondary N-factors must address the issues of mode selection from the competing instability 
mechanisms as well as the nonlinear interactions between multiple types of primary and secondary instabilities and 
the breakdown scenarios. An additional key question in the context of transition prediction for laminar flow 
technology corresponds to whether or not increasingly rational transition prediction can be achieved with easily 
available input data regarding the disturbance environment such as the surface finish and a small number of 
parameters characterizing the freestream disturbances. 
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